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NPSNET:  Flight  Simulation  Dynamic 
Modeling  Using  Quaternions 


Abstract 

The  Naval  Postgraduate  School  (NPS)  has  actively  explored  the  design  and  implementa¬ 
tion  of  networked,  real  time,  three-dimensional  battlefield  simulations  on  low-cost, 
commercially  available  graphics  workstations.  The  most  recent  system,  NPSNET,  has 
improved  in  functionality  to  such  an  extent  that  it  is  considered  a  low-cost  version  of  the 
Defense  Advanced  Research  Project  Agency’s  (DARPA)  SIMNET  system.  To  reach  that 
level,  it  was  necessary  to  economize  in  certain  areas  of  the  code  so  that  real  time  per¬ 
formance  occurred  at  an  acceptable  level.  One  of  those  areas  was  in  aircraft  dynamics. 
However,  with  “off-the-shelf”  computers  becoming  faster  and  cheaper,  real-time  and 
realistic  dynamics  are  no  longer  an  expensive  option.  Realistic  behavior  can  now  be 
enhanced  through  the  incorporation  of  an  aerodynamic  model.  To  accomplish  this  task, 
a  prototype  flight  simulator  was  built  that  is  capable  of  simulating  numerous  types  of 
aircraft  simultaneously  within  a  virtual  world.  Besides  being  easily  incorporated  into 
NPSNET,  such  a  simulator  also  provides  the  base  functionality  for  the  creation  of  a  gen¬ 
eral  purpose  aerodynamic  simulator  that  is  particularly  useful  to  aerodynamics  students 
for  graphically  analyzing  differing  aircraft’s  stability  and  control  characteristics.  This  sys¬ 
tem  is  designed  for  use  on  a  Silicon  Graphics  workstation  and  uses  the  GL  libraries.  A 
key  feature  of  the  simulator  is  the  use  of  quaternions  for  aircraft  orientation  representa¬ 
tion  to  avoid  singularities  and  high  data  rates  associated  with  the  more  common  Euler 
angle  representation  of  orientation. 


I  Introduction 

The  current  state  of  the  art  in  simulation  technology  has  provided  today’s 
military  with  many  valuable  training  experiences  that  could  not  have  been  ob¬ 
tained  elsewhere  and,  as  a  result,  has  gready  increased  survivability  and  readi¬ 
ness.  From  flight  simulators,  which  allow  a  pilot  to  explore  the  edge  of  the 
flight  envelope  without  endangering  crew  or  multimillion  dollar  assets,  to  bat¬ 
tlefield  simulators,  which  allow  entire  fighting  divisions  to  practice  command 
and  control  without  having  to  incur  the  enormous  costs  of  running  a  full  blown 
field  exercise,  computer  simulation  has  become  a  way  of  doing  business  within 
the  military. 

One  simulation  system  designed  by  the  Defense  Advanced  Research  Projects 
Agency  (DARPA)  is  the  Simulation  Networking  system  (SIMNET)  (Thorpe, 
1987).  SIMNET  is  a  networked  battlefield  simulator  that  allows  multiple  user 
interaction  on  the  battlefield  at  many  different  levels.  Vehicle  simulators,  such 
as  tanks  and  aircraft,  connect  to  the  network  and  become  part  of  a  three-dimen¬ 
sional  world.  At  the  Naval  Postgraduate  School  (NPS),  an  effort  to  develop  a 
SIMNET-type  system  based  on  commercially  available,  general  purpose, 
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graphic  workstations  has  been  active  for  a  number  of 
years.  This  system,  NPSNET,  consists  of  Silicon  Graph¬ 
ics  workstations  attached  to  a  local  area  Ethernet  (Zyda, 
Pratt,  Monahan,  &  Wilson,  1992).  Eventually, 
NPSNET  will  become  a  node  on  the  SIMNET  network. 

The  speed  of  the  computer  platforms  on  which 
NPSNET  now  runs  has  increased  significantly  since  its 
inception.  It  soon  became  evident  that  more  realistic 
vehicle  dynamics  was  desirable  and  would  substantially 
improve  system  behavior.  Therefore,  this  work  is  the 
result  of  the  research  done  and  the  methodology  used 
for  providing  this  additional  functionality  to  NPSNET’s 
aircraft  simulations. 

There  are  several  issues  to  be  addressed  when  incorpo¬ 
rating  an  aerodynamic  model  into  a  computer  simula¬ 
tion.  The  complexity  of  the  aerodynamic  model,  which 
orientation  model  to  use,  and  how  aircraft  data  should 
be  represented  in  the  system  are  the  critical  issues  and 
are  the  center  of  focus  in  this  work.  Of  primary  impor¬ 
tance  is  that  the  complexity  of  the  model  fits  the  objec¬ 
tive  of  the  simulation.  A  complete  aerodynamic  model 
that  includes  fully  articulated  control  surfaces  and  air¬ 
flow  divergence  patterns  over  the  aircraft  would  seri¬ 
ously  affect  real  time  performance  on  any  computer. 
Such  models  are  usually  computed  in  non-real  time  on 
supercomputers  and  are  not  appropriate  for  use  on  low- 
cost  graphics  workstations.  On  the  other  hand,  model¬ 
ing  the  dynamics  of  an  aircraft  kinematically,  so  that  the 
aircraft’s  velocity  and  orientation  are  a  linear  result  of 
control  input,  does  not  reproduce  the  nuances  of  aircraft 
motion  and  response  that  a  user  of  a  flight  simulation 
would  expect.  The  aerodynamic  model’s  complexity 
must  provide  as  much  realism  as  possible  without  reduc¬ 
ing  the  frame  rate  below  an  acceptable  level. 

The  choice  of  orientation  model  considered  is  be¬ 
tween  the  Euler  angle  or  quaternion  approach.  Which 
one  to  use  has  been  the  subject  of  heated  debate  among 
computer  scientists  (Goldiez  &  Lin,  1991).  Either 
model  can  be  used  to  represent  orientation,  but,  depend¬ 
ing  on  the  simulation’s  objectives,  one  method  has  cer¬ 
tain  advantages  over  the  other.  It  basically  comes  down 
to  determining  which  approach  provides  the  right  tool 
for  the  job  (Shoemake,  1985). 

Because  each  type  of  aircraft  exhibits  its  own  specific 


aerodynamics  and  handling  characteristics,  it  is  desirable 
to  change  these  characteristics  depending  on  the  type  of 
aircraft  simulated.  One  solution  is  to  base  the  aerody¬ 
namic  model  on  the  aircraft’s  stability  coefficients,  iner¬ 
tial  coefficients,  and  airframe  specifications,  all  of  which 
are  available  in  most  aerodynamic  stability  and  control 
textbooks.  Stability  coefficients  provide  a  very  accurate 
model  of  aircraft  flight  behavior.  However,  care  must  be 
taken  when  modeling  some  of  the  newer  generation 
fighters.  To  improve  maneuverability,  these  aircraft  have 
been  designed  aerodynamically  neutral  to  unstable. 

Their  stability  coefficients  reflect  this  instability. 

2  Coordinate  Systems  and  Terminology 

Coordinate  systems  and  the  method  in  which  they 
are  described  vary  to  a  great  extent  depending  on  the 
application  and  the  preference  of  the  user.  In  aircraft 
simulations,  coordinate  systems  fall  into  two  broad 
classes,  “body”  coordinates  and  “earth”  or  “inertial”  co¬ 
ordinates  (Rolfe  &  Staples,  1986).  Body  coordinates 
have  their  origin  based  at  the  center  of  gravity  and  con¬ 
tinually  move  with  the  aircraft.  Inertial  coordinates,  on 
the  other  hand,  are  defined  with  respect  to  the  earth  and 
have  their  origin  positioned  at  some  suitable  location 
such  as  the  center  of  the  simulated  world.  Other  coordi¬ 
nate  systems  exist,  based  on  parameters  such  as  the  flight 
path  and  angle  of  attack.  However,  the  aerodynamic 
model  presented  in  this  paper  is  based  on  geometric 
body  and  world  coordinate  systems.  In  general,  all  aero¬ 
dynamic  forces,  accelerations,  and  velocities  are  calcu¬ 
lated  in  the  body  coordinate  system  first,  and  then  con¬ 
verted  to  the  world  coordinate  system  prior  to  updating 
an  aircraft’s  position  and  altitude. 

Figure  1  shows  the  generally  accepted  convention  for 
labeling  of  the  axes  in  the  two  coordinate  systems 
(Anderson,  1989).  Body  coordinates  are  defined  with 
the  origin  at  the  center  of  gravity  (CG),  the  x  axis  along 
the  fuselage  pointing  out  the  nose  of  the  aircraft,  the  y 
axis  along  the  wing-line  pointing  out  the  right  wing,  and 
the  z  axis  pointing  out  the  bottom  of  the  plane.  World 
coordinates  are  defined  with  the  origin  based  at  a  fixed 
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Body  Coordinates 


World  Coordinates 


Zw 


[w 


North 


■^East 


Down 


xw 


Figure  I .  World  and  body  coordinate  systems. 


Xw. 


Y  7 

1  W’ 


Aircraft  location  in  world  coordinates  (feet) 


Uw,  Vw,  Ww  Aircraft  velocity  in  world  coordinates  (ft/sec) 

\| /,  0,  0  Azimuth,  Elevation,  Roll  in  world  coordinates  (radians) 


Figure  2.  Terms  defined  within  the  world  coordinate  system. 


U,  V,W 
P,  Q,R 


U,  V,  w 
P.Q.R 

Fx>  Fy-  Fz 

L,M,N 

a 

P 


Linear  velocity  along  X,  Y,  and  Z  body  axes  (ft/sec) 
Angular  velocity  around  X,Y,  and  Z  body  axes  (rad/sec) 

Resultant  velocity  Vector  ^2  +  y2  +  w2 

Wind  velocity  across  tail  of  aircraft 

Linear  acceleration  (ft/sec2) 

Angular  acceleration  (rad/sec2) 

Forces  acting  on  aircraft 
Moments  about  the  X,  Y,  and  Z  axes 
Angle  of  attack[tan_1  (W/U)] 

Sideslip  [tan'1  (V/U)] 


Figure  3.  Terms  defined  within  the  aircraft  body  coordinate  system. 


point  on  the  ground,  thtx  axis  pointing  north,  they  axis 
pointing  east,  and  the  z  axis  pointing  down.  Because  of 
its  limited  effect,  the  curvature  of  the  earth  is  usually  ig¬ 
nored. 

In  a  dynamics  model,  velocities,  accelerations,  and 
forces  are  described  in  both  world  and  body  coordinate 
systems.  Without  an  explicit  description  of  the  variables 
used  to  describe  the  model,  confusion  can  arise.  Most 
terms  described  in  this  paper  refer  to  the  geometric  body 
axes.  However  if  a  reference  is  made  to  the  world  coor¬ 
dinate  system  the  subscript  tcw”  is  used  (Fig.  2). 

Terms  without  the  ccw”  subscript  relate  to  body  axes 
and  include  linear  and  angular  velocities,  accelerations, 


forces,  moments,  angle  of  attack,  and  sideslip  angle  (Fig. 
3).  Note  that  the  direction  of  angular  accelerations  and 
velocities  and  moment  terms  are  defined  using  the  right- 
hand  rule  around  their  respective  axes  (Fig.  4). 

The  aircraft  control  surfaces  such  as  elevator,  ailerons, 
and  rudder  are  defined  as  a  rotation  in  radians  around 
their  respective  hinge  points  on  the  aircraft.  When  a  con¬ 
trol  surface  is  flush  with  the  aircraft,  the  angle  of  deflec¬ 
tion  is  zero  (Fig.  5). 

Within  the  aerodynamic  model,  the  particular  aircraft 
being  modeled  is  characterized  by  certain  dimensional 
characteristics.  A  description  of  these  terms  is  included 
in  Figure  6. 


8e 


8a 

5r 


Figure  5. 


Figure  4.  Notation  with  respect  to  body  axes. 


3  Aerodynamic  Model 

The  mathematical  model  presented  takes  forces, 
control  inputs,  and  aircraft  specifications  as  inputs,  gen¬ 
erating  linear  and  angular  velocities  in  aircraft  body  co¬ 
ordinates  as  outputs  (Fig.  7).  Based  on  a  classical  repre¬ 
sentation  of  linear  aerodynamics  and  utilizing  the  total 
force  Eqs.  (3.18)-(3.20),  all  forces  associated  with  lift 
and  drag  are  calculated  utilizing  aerodynamic  stability 
derivatives  (Roskam,  1979). 

Stability  derivatives,  first  used  over  a  half-century  ago, 
assume  that  all  aerodynamic  forces  and  moments  can  be 
expressed  as  a  function  of  the  instantaneous  value  of  the 
perturbation  variables  (Nelson,  1989).  The  perturbation 
variables  are  the  instantaneous  changes  from  the  refer¬ 
ence  conditions  of  translational  velocities,  angular  veloci¬ 
ties,  control  deflections,  and  their  derivatives.  For  exam¬ 
ple,  the  term  8X/8 u  is  the  stability  derivative  defining 
the  change  inX  force  with  respect  to  the  change  in  for¬ 
ward  speed.  This  derivative  can  be  expressed  in  terms  of 
a  nondimensional  coefficient  CXu  as  follows: 

8X  1 

u  =  c*-«fis  <3'‘> 

where 


elevator  deflection  positive  down  (radians). 

A  positive  8e  produces  a  positive  lift  and  a  negative  pitch  moment. 

aileron  deflection  positive  left  (radians). 

A  positive  8a  produces  a  negative  roll  moment. 

positive  nose  left  (radians).A  positive  8r  produces 
a  positive  sideforce  and  a  negative  yaw  moment. 

Terminology  defining  aircraft  controls. 


Figure  8  lists  the  nondimensional  coefficients  used  in 
this  model.  These  coefficients  are  generally  broken  down 
into  three  categories,  lateral,  longitudinal,  and  control. 
The  longitudinal  coefficients  represent  forces  effecting 
the  longitudinal  axes  of  the  aircraft,  while  the  lateral  co¬ 
efficients  represent  forces  affecting  the  lateral  axes  of  the 
aircraft.  Nondimensional  coefficients,  generated  in  actual 
aircraft  testing,  are  available  for  most  aircraft.  By  using 
these  coefficients  in  combination  with  the  dynamics 
equations,  it  is  possible  to  build  a  general  use  flight  sim¬ 
ulator. 

Using  the  nondimensional  coefficients,  lift,  drag,  and 
sideforce  are  calculated  as  follows: 


V  = 


Is  Is 

Clo  CjLa(x  ^  C^cl  ^ 

(Ft  +  AFe)l2 


D  = 


+  Q.seSe 


CIX)  +  CDa  a  +  CDil,be 


Ft 

(Ft  +  AFe)l2 


P  Ft2S 

—  <3-3) 


Ft 


P  v2ts 


S¥  -  [Crp(3  +  CysrSr] 


P  V27S 


(3.4) 

(3.5) 


Once  lift,  drag,  and  sideforce  are  calculated,  these 
forces  are  translated  into  forces  along  the  aircraft  X,  T, 
andZ  axes  as  shown  in  Eqs.  (3.6)  through  (3.8).  The 
terms  Ea*,  EAr,  and  F^z  represent  the  resultant  aerody¬ 
namic  forces. 


SC* 

Cxu  ~  8(«/«0)  (3‘2 

and  Q  is  the  dynamic  pressure,  VipU?,  where  p  is  the  air 
density  at  the  aircraft  altitude. 


Faz  “  -L '  cos  a  -  D  sin  a  (3.6) 

Fax  -  L '  sin  a  -  D  cos  a  -  SF  sin  p  (3.7) 

FAy  =  cos  p  (3.8) 
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S  surface  area  of  wing  (ft2) 

b  wing  span  (ft) 

c  chord  length  (ft) 

w  weight  (lbs) 

Ixx  roll  inertia  (slug-ft2) 

Iyy  pitch  inertia  (slug-ft2) 

Izz  yaw  inertia  (slug-ft2) 

Figure  6.  Aircraft  dimensional  specifications. 


Figure  7.  Basic  aerodynamic  model. 


The  aerodynamic  moments  represent  the  torque 
forces  about  the  center  of  the  aircraft  and  are  determined 
in  the  following  equations: 


La  = 


Qp(3  +  CLPP  2Vt  +  ClrR  2yT 


+  Cx5a8a  +  C/5,f)r 


pVt  2Sb 


(3.9) 


Ma  = 


i/ 

4  T  Cmq(^  ^Vt 


4  Cm6l 2J/j  C^cSe 


(Ft  4  AFe) 


Ft 


pVr2Sc 


(3.10) 


Na  = 


b  b 

4  CNPP  4  CnrR  2^ 


+  C^5a8a  +  Cn  5r8r 


pFr25& 

—  <3-n> 


The  forces  and  moments  that  result  from  the  above 
calculations  are  added  to  other  forces  and  moments  at 
this  time: 


Lx  —  Lax  +  ^Thrust 

(3.12) 

Lr  =  Lat 

(3.13) 

Lz  =  LAz 

(3.14) 

L  La  4  Torque 

(3.15) 

JV  —  LA.a  4  I^L-Thrust  4  iVf(3yro 

(3.16) 

H  Ha  4  ^Thrust  4  NQyro 

(3.17) 

Engine  forces  such  as  thrust,  torque,  and  gyroscopic 
effect  as  well  as  environmental  forces  such  as  wind  shear 
can  have  anywhere  from  a  minor  to  a  significant  effect 
on  the  forces  and  moments  along  all  axes  of  the  aircraft 
(Roskam,  1979).  However,  to  limit  this  complexity  of 
the  model,  some  simplifications  are  made.  Engine  thrust 
is  limited  to  the  X-axis  only  and  no  calculations  are  made 
for  torque  or  gyroscopic  effect  since  one  of  the  author’s 
experiences  as  a  pilot,  and  reference  to  the  relevant  litera¬ 
ture,  indicates  that  these  are  second-order  effects  for 
high-performance  aircraft. 

The  total  force  equations  are  used  to  determine  the 
linear  acceleration  of  the  aircraft  (Nelson,  1989): 

Lx 

U  =  VR  -  WQ  -g  sin  0  4  ^  (3.18) 

Fr 

V  =  WP  -  UR  Fa  sin  cj>  cos  0  -I - 

^  m 

(3.19) 
Fz 

W  =  UQ  -  VP  Fa  cos  cj)  cos  0  H - 

^  ^  m 

(3.20) 

The  total  moment  equations  are  used  to  derive  the  equa¬ 
tions  for  solving  angular  acceleration: 

L  =  IxxP  -  IxzR  -  IxzPQ  +  (Izz  -  Irr)RQ  (3.21) 

M  =  IrrQ  +  (Ixx  ~  hz)PR  +  {hz{P2  ~  R2)  (3-22) 

N  =  IZZR  -  IXZP  +  (Irr  -  I^PQ  +  IXZQR  (3.23) 

However,  prior  to  solving  for  either  P  or  R,  an  interim 
step  is  required: 

L"  =  L  +  IXZPQ  -  (Izz  -  Irr)RQ  (3.24) 


Longitudinal  Coefficients 


"Lo 


"Do 


"La 


"Da 


"Mo 


"Ma 


LQ 

"MQ 


La 

C  . 

Ma 

Lateral  Coefficients 

CYp 

CLP 

CLP 

CLR 

CN(3 

CNP 

CNR 

Control  Coefficients 


C 

C 

C 

C 


L8e 

D8e 

M8e 

L8a 


Reference  Lift  at  zero  angle  of  attack 
Reference  Drag  at  zero  angle  of  attack 
Lift  curve  slope 
Drag  curve  slope 
Pitch  moment 

Pitch  moment  due  to  angle  of  attack 

Lift  due  to  pitch  rate 

Pitch  moment  due  to  pitch  rate 

Lift  due  to  angle  of  attack  rate 

Pitch  moment  due  to  angle  of  attack  rate 

Side  force  due  to  sideslip 
Dihedral  effect 
Roll  damping 
Roll  due  to  yaw  rate 
Weather  cocking  stability 
Rudder  adverse  yaw 
Yaw  damping 

Lift  due  to  elevator 
Drag  due  to  elevator 
Pitch  due  to  elevator 
Roll  due  to  aileron 


Figure  8.  Aircraft  specification  notation. 


N'=N-(Irr-  I^PQ  -  IXZRQ 

(3.25) 

Therefore, 

P  =  (L "Izz  -  N'IXZ)/(IXXIZZ  -  I2XZ) 

(3.26) 

Q  =  (M  —  (- 1 xx  ~  Izz)PR 

-  IXZ(P 2  -  R2))Hrr 

(3.27) 

R  =  (N'lxx  +  L"IXZ)/(IXXIZZ  -  I\z) 

(3.28) 

are  the  equations  for  angular  acceleration. 

Linear  and  angular  velocities  are  determined  by  nu¬ 
merically  integrating  the  accelerations.  The  trapezoidal 
rule,  sometimes  referred  to  as  the  modified  Euler 
method  or  the  first-order  predictor-corrector  method,  is 


used.  The  general  method  for  this  integration  technique 
is  as  follows  (Press,  Flannery,  Teukolsky,  &  Vetterling, 
1990): 

(Pn-i  +  Pn\ 

Pn=Pn- 1  +  [ - 2 - )dt  (3'29) 

where 

Pn  —  new  value  of  P 
Pn~i  =  previous  value  ofP 
Pn  =  predicted  rate-of-change  of  P 
Pn~\  =  previous  rate-of-change  ofP 
dt  =  integration  step  size 
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4  Position  Updates 


Because  position  updates  usually  occur  in  a  world 
coordinate  system,  the  aircraft’s  linear  velocities  must  be 
converted  into  world  position  rates  by  applying  the  fol¬ 
lowing  transformation,  which  effectively  rotates  the 
(U  V  W)  vector  by  the  Euler  angles  (Roskam,  1979). 
The  order  of  transformation  is  important  when  utilizing 
this  method  and  proceeds  as  follows: 


U, 

1  0 

0 

u 

n 

— 

0  COS 

<t>  ■ 

-sin  <j> 

V 

0  sin 

COS  (|) 

\ w 

cos  0 

0 

sin  0 

A' 

V* 

= 

0 

l 

0 

n 

W»_ 
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Integrating  the  resultant  velocity  vector,  now  in  world 
coordinates,  by  the  time  step  of  the  program,  a  position 
update  is  obtained  [Eq.  (4.4)]: 
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This  technique  is  revisited  in  the  next  section. 


5.  Methods  of  Orientation 

The  aerodynamic  model  generates  rotational  veloc¬ 
ities  relative  to  the  fixed  aircraft  body  coordinate  system. 
But  just  as  position  updates  could  only  be  determined 
after  converting  linear  velocities  into  world  coordinates, 
orientation  updates  require  conversion  of  angular  veloci¬ 
ties  in  a  similar  manner.  Three  methods  exist  for  defin¬ 
ing  the  conversion  of  angular  velocities  to  orientations  in 
world  coordinates,  each  having  its  own  particular  advan¬ 
tages  and  disadvantages  (Goldiez  &  Lin,  1991). 


The  most  popular  of  these  three  methods  is  known  as 
the  Euler  method.  Using  a  sequence  of  three  angles,  the 
Euler  method  provides  an  intuitive  description  of  air¬ 
craft  attitude  in  world  space  (Rolfe,  1986).  These  angles 
consist  of  the  familiar  azimuth  angle  i| the  elevation 
angle  0,  and  the  roll  angle  <J>.  The  next  method,  which 
has  become  popular  in  recent  years,  is  the  quaternion 
method.  Based  on  the  unit  sphere,  the  quaternion 
method  provides  an  elegant  method  of  defining  rota¬ 
tions  through  the  use  of  four  parameters.  Three  of  the 
coordinates  describe  the  axis  of  rotation  while  the  fourth 
is  determined  by  the  angle  through  which  the  rotation 
occurs  (Shoemake,  1985).  The  third  method  of  defining 
orientation  is  the  direction  cosine  matrix.  The  direction 
cosines  relate  the  aircraft  body  axis  frame  to  the  world 
reference  frame.  Direction  cosines,  as  used  in  flight  sim¬ 
ulation,  are  generally  determined  from  either  Euler  an¬ 
gles  or  quaternions  and  are  utilized  for  transformations 
between  axes.  However,  an  alternative  approach  (not 
discussed  in  this  paper)  is  to  use  incremental  rotation 
matrices  to  update  rotation  matrices  (Paul,  1981).  A 
disadvantage  of  this  approach  is  that  repeated  incremen¬ 
tal  rotation  matrix  multiplication  can  result  in  drift  re¬ 
quiring  periodic  renormalization  of  the  direction  cosine 
matrix  (Funda,  Taylor,  &  Paul,  1990). 

Each  method  has  its  own  particular  advantages  and 
disadvantages  and  their  use  depends  on  the  application 
and  the  implementation.  Because  NPSNET  is  a  net¬ 
worked  simulator,  the  orientation  model  used  must  not 
only  render  orientations  in  the  world  of  the  aircraft  be¬ 
ing  piloted,  but  also  of  other  aircraft  in  the  world,  either 
flying  autonomously  or  piloted  remotely  across  the  net¬ 
work. 


6  Euler  Method 

The  most  common  method  of  defining  an  aircraft’s 
orientation  in  world  space  is  by  the  Euler  attitude  angles. 
Starting  with  the  aircraft’s  axis  origin  aligned  with  the 
world’s  axis  origin,  the  Euler  angles  specify  three  succes¬ 
sive  rotations  to  bring  the  world  coordinates  into  align¬ 
ment  with  the  aircraft.  The  fact  that  there  exist  12  possi¬ 
ble  ways  to  define  rotations,  each  with  potentially 


rotate  in  \j /  around  Z 
Z 


rotate  in  0  around  Y 


Figure  9.  Euler  attitude  angle  rotation. 


different  results,  means  that  the  order  of  these  rotations 
is  important.  Euler  chose  the  convention  of  rotating  first 
about  the  z  axis,  then  about  the  new#  axis,  and  finally 
about  the  new  z  axis.  This  convention  exists  in  celestial 
mechanics,  applied  mechanics,  and  molecular  and  solid 
state  physics.  The  convention  used  in  quantum  mechan¬ 
ics,  nuclear  physics,  and  particle  physics,  chooses  to  ro¬ 
tate  first  about  the  z,  then  the  newjy,  and  finally  the  new 
z  (Burchfiel,  1990).  The  convention  most  often  used  in 
graphics  is  standard  to  aerospace  engineers  and  has  been 
proposed  for  use  by  SIMNET  (UCF/IST  1990).  Using 
the  right-hand  rule,  rotations  are  made,  first,  about  the  z 
axis  by  the  angle  i| /,  then  about  the  newy  axis  by  angle  0, 
and  finally  about  the  new#  axis  by  angle  c|>  (Fig.  9). 

The  range  of  values  the  attitude  angles  can  take  are 

IT 

l|i  =  ±7T  0  =  ±—  (j)  =  ±TT 

The  aerodynamic  model  generates  velocities  in  body 
coordinates.  As  seen  in  Eqs.  (4.1)-(4.3),  the  linear  body 
rates  are  transformed  into  world  rates  by  application  of 
the  Euler  angles.  What  follows  is  a  method  for  obtaining 
these  angles. 

There  is  a  direct  relationship  between  Euler  attitude 
angles  and  the  angular  velocity  of  the  aircraft  around  its 
body  axes  (Nelson,  1989).  From  this  relationship,  the 
rates  of  change  of  the  attitude  angles  can  be  derived: 

cj>  =  P  +  Q  sin  <\>  tan  0  +  R  cos  $  tan  0  (6.1) 

0  =  Q  cos  <\>  -  R  sin  <\)  (6.2) 

i|i  =  Q  sin  <\>  sec  0  +  R  cos  <\>  sec  0  (6.3) 

The  inverse  of  the  above  equations  are 

(6.4) 


rotate  in  <|>  around  X 


<£)  =  0  cos  4>  -t-  vj/  sin  4>  cos  0  (6.5) 

R  =  —0  sin  +  iji  cos  4>  cos  0  (6.6) 

Equations  (6.1)  through  (6.3)  are  also  known  as  the 
gimbal  equations  and  are  quite  commonly  used  in  simu¬ 
lation.  However,  a  problem  exists  when  pitch,  0,  goes 
through  the  vertical.  That  is,  where  pitch  becomes 
±  (tt/2).  At  that  point  c)>  and  i|;  become  undefined.  Im¬ 
plementing  a  flight  dynamics  model  capable  of  complete 
vertical  maneuvering  necessitates  “fixing”  the  code  so  a 
division  by  zero  does  not  occur. 


7  Direction  Cosines 


In  the  case  of  a  flight  simulation,  transforming  be¬ 
tween  body  coordinates  and  world  coordinates  is  done 
quite  frequently.  A  convenient  way  to  represent  the 
transformation  between  two  coordinate  systems  is  with 
the  direction  cosines.  Using  matrix  notation  and  the  di¬ 
rection  cosines  (a,  b ,  c),  the  transformation  from  body  to 
world  axes  is  expressed  by 
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X,  T,  and  Z  represent  vectors  of  any  kind,  such  as 
force,  velocity,  and  acceleration.  The  inverse  relation¬ 
ship,  converting  world  coordinates  to  body  coordinates, 
is  the  transpose: 
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=  cos0cos\|/  bj  =  sin(|)sin0cos\|/- cos(|)sin\|/ 

a2  =  cos0sin\}/  b2  =  sin^sinGsinv;/  +  cos(J)cos\|/  Cr 

a^  =  —  sin0  =  sin<|)cos0 

Figure  1 0.  Direction  cosines  in  terms  of  Euler  angles. 

In  terms  of  the  Euler  attitude  angles,  the  direction 
cosines  for  the  above  transformations  are  shown  in  Fig¬ 
ure  10. 

It  should  be  noted  that  there  are  12  ways  in  which 
Euler  angles  can  be  defined,  and,  as  a  result,  just  as  many 
ways  to  compute  the  direction  cosines,  although  the  val¬ 
ues  finally  obtained  are  independent  of  the  choice  of  Eu¬ 
ler  angles. 

The  direction  cosines  are  needed  for  transformations 
between  coordinate  systems,  whether  Euler  angles  or 
quaternions  are  used  to  specify  orientation.  [Although, 
as  explained  in  Funda  et  al.  (1990),  an  alternative  to  ma¬ 
trix  multiplication  using  direction  cosines  is  to  transform 
points  by  quaternion  multiplication.]  Direction  cosines 
were  already  used  for  transforming  linear  velocities  in 
body  coordinates  to  world  coordinates.  By  multiplying 
the  transformation  matrices  of  [Eqs.  (4. 1)— (4.3)]  into 
one  matrix,  the  result  would  be  identical  to  the  transfor¬ 
mation  matrix  of  Eq.  (7.1).  By  using  direction  cosines, 
the  need  for  determining  the  intermediate  velocities  is 
eliminated. 

8  Quaternion  Method 

An  alternate  method  that  has  gained  popularity  in 
the  graphics  community  in  the  mid-1980s  is  through  the 
use  of  the  unit  quaternion.  Not  a  new  method,  quaterni¬ 
ons  have  been  around  for  over  a  century.  Augmenting 
the  “four-parameter  method,”  they  have  been  useful  to 
aerodynamic  engineers  for  some  time  and  are  still  the 
method  of  choice  for  describing  spacecraft  orientation 
(Mitchell  Sc  Rodgers,  1965).  Discovered  by  Sir  William 
Rowan  Hamilton  in  1843  as  a  result  of  a  search  for  a 
generalization  of  complex  numbers,  quaternions  provide 
an  efficient  means  for  updating  orientations  (Shoemake, 
1985). 

There  are  numerous  ways  to  interpret  the  quaternion 
mathematically.  They  can  be  described  as  an  algebraic 


=  cos<|)sin0cos\j/4-  sin<|)sin\|r 
=  cos<()sin0sin\|/-  sin([)cos\|/ 
Cg  =  cos  (j)  cos  0 


Figure  1 1 .  Quaternion  orientation. 

quantity, 

w  +  ix  +  jy  +  kz  (8-1) 

as  a  point  in  three-dimensional  projective  space 
(w>,  x,y,  2),  as  a  linear  transformation  of  four  space  (ma¬ 
trix),  or  as  a  scalar  plus  3 -vector: 

( w ,  v)  v  =  ix  +  jy  +  kz  (8.2) 

The  best  notation  depends  on  their  intended  use.  The 
most  intuitive  approach  is  to  view  the  quaternion  as  a 
scalar  plus  3-vector  [Eq.  (8.2)].  However,  for  algebraic 
manipulation,  Eq.  (8.1)  generally  becomes  more  useful. 

A  common  way  of  defining  quaternion  orientation  is 
in  combination  with  Euler’s  theorem  which  states  that 
the  orientation  of  a  rigid  body  can  be  described  as  a  rota¬ 
tion  about  an  axis  v  by  rotation  angle  (Fig.  11)  (Gold¬ 
stein,  1980).  Constraining  the  axis  vector  v  to  be  of  unit 
magnitude,  the  quaternion  becomes 

<E  O 

Q  =  cos  —  ’  v  sin  —  (8.3) 

This  representation  is  always  of  unit  magnitude  such 
that 

w2  +  x2  +y2  +  z2  =  1  (8.4) 

Prior  to  defining  how  to  rotate  a  rigid  body  using  the 
unit  quaternion,  it  is  necessary  to  review  some  of  the 
mathematics  associated  with  the  quaternion.  For  the 
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Figure  12.  Four  parameter  method. 


purpose  of  a  flight  simulation,  an  understanding  of  the 
multiplication  and  the  quaternion  derivative  is  necessary. 
More  comprehensive  reviews  are  available  in  Shoemake 
(1985),  Goldstein  (1980),  Funda  et  al.  (1991),  and 
Chou  (1992). 

Most  applications  involving  quaternions  make  use  of 
the  mathematics  associated  with  their  multiplication. 
Similar  to  the  algebra  associated  with  imaginary  num¬ 
bers,  quaternions  have  three  imaginary  units,  ij\  and  k 
and  are  noncommunitive  under  multiplication  with 

i 2  =j2  =  k2=  -1  (8.5) 

and 

ij  =  k  =  —ji  jk  =  i  =  —kj  ki  —  j  —  —ik  (8.6) 

In  algebraic  notation,  the  product  of  quaternion  Q  mul¬ 
tiplied  by  quaternion  is 

QQ\  =  (W  +  ix  +jy  +  kz)(w1  +  ixy  +  jyx  +  kzx) 

=  (ww\  -  xxl  —  yj\  -  zzY) 

+  i(xw i  +  wx1  -  zyi  +  yzi)  (8.7) 


The  results  of  the  above  multiplication  is  a  rotation 
from  the  orientation  represented  by  Q  to  the  new  cumu¬ 
lative  orientation  of  Q  and  in  quaternion  terms. 

Multiplication  provides  a  method  of  orientation  ex¬ 
trapolation  that  can  be  of  benefit  in  a  networked  simula¬ 
tion.  IfjQi  represents  a  finite  rotation  based  on  an  inte¬ 
gral  time  step,  and  Q  represents  the  cumulative  rotation, 
then  a  repeated  multiplying  of  Q  by  Qi  will  result  in  a 
smooth  rotation  across  a  series  of  update  frames  (Burch- 
fiel,  1990). 

One  frame  of  axes  (body  coordinates)  can  be  brought 
into  coincidence  with  a  reference  frame  by  a  single  rota¬ 
tion  D  about  a  fixed  axis  making  angles  ^4,  B,  and  C  with 
a  second  reference  frame  (world  coordinates).  The  four 
parameters^!,  B,  C,  and D,  therefore,  define  the  orienta¬ 
tion  of  the  aircraft  body  in  world  coordinates  (Rolfe  & 
Staples,  1986).  The  transformation  matrix  relating  body 
to  world  coordinates  using  these  four  parameters  is 
shown  in  Figure  12.  While  this  matrix  involves  four  an¬ 
gles  and  appears  to  be  more  complex  than  the  Euler  an¬ 
gle  matrix  of  Figure  10,  it  can  be  simplified  by  making 
the  substitutions: 


+  j{yv?\  +  zx i  -  wyi  +  xzi) 
+  k(zwi  +  yxi  —  xyi  +  wz{) 
In  vector  notation: 


q0  =  cos  1/2  D 
qY  =  cos  A  sin  1/2D 
q2  =  cos  B  sin  1/2 D 


(8.9) 


QQ\  =  (W,  v)(wu  yl)  =  wwl-v-  vl5 


wv1  +  w\\  +  v  x  vx  (8.8) 


qz  =  cos  C  sin  1/2JD 


The  transformation  matrix  then  becomes 
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which  represents  the  transform  based  on  the  unit  quater¬ 
nion.  This  result  is  used  for  determining  position  up¬ 
dates  as  well  as  orientation  updates. 

In  case  the  angles^,  B ,  C,  D  are  not  initially  known, 
values  for  q2 ,  #3,  and  q4  can  be  obtained  by  a  straight¬ 
forward  process  from  the  initial  direction  cosine  matrix 
(Funda  et  al.,  1991).  To  update  the  resulting  quaternion 
from  angular  accelerations,  the  following  equations  are 
used: 


q0  =  -l/2(qiP  +  q2Q  +  q3R) 
9i  =  l/2(0o P  +  92R  ~  93Q) 
q2=  l/2(q0Q  +  qzP  -  qxR) 
qz  =  1/2 (q0R  +  qxQ-  q2P) 


Because  of  the  constraint  that  the  quaternion  be  of  unit 
value  and  assuming  an  integration  step  size  of  less  than 
1,  the  above  set  of  equations  become 


q0  =  -1/2  (qYP  +  q2Q  +  q3R)  +  \q0 

h\  =  l/2(9op  +  92 R  ~  <hQ)  +  Mi 

q2  =  1/2(^0j2  +  <hp  “  9ir)  +  M2 
q3  =  l/2(q0 R  +  qxQ  -  q2P)  +  \q3 


(8.12) 


where  \  is  an  integration  drift  correction  gain  given  by 


^  =  1  _  (?o  +  9i  +  9\  +  ?i)  (8.13) 


Alternatively,  Eq.  (8.11)  can  be  integrated  without  drift 
correction  providing  that  periodic  normalization  to  unit 
magnitude  is  accomplished  (Funda  et  al.,  1990). 

Many  of  the  auxiliary  computations  involved  with  a 
flight  simulation  require  the  use  of  Euler  angles.  It  is  to 
be  emphasized,  however,  that  knowledge  of  Euler  angles 
is  not  required  to  obtain  the  direction  cosines  of  Eq. 
(8.10).  In  fact,  Eq.  (8.11)  can  be  used  to  update  the  di¬ 


rection  cosine  matrix  without  the  calculation  of  any  trig¬ 
onometric  functions  at  all. 

When  needed,  Euler  angles  can  be  obtained  from  the 
transformation  matrix  in  Eq.  (8.10)  via  the  following 
method.  Because  pitch  is  limited  to  ±tt/2,  cos  (0)  is  al¬ 
ways  positive.  As  a  result,  obtaining  these  angles  is  rela¬ 
tively  simple.  The  evaluation  angle  is  derived  from  the 
transformation  matrix  [Eq.  (7.1)]  as  follows: 

0  =  asin  (—03)  (8.14) 

To  obtain  the  azimuth  (1 |i)  angle  it  must  be  noted  that, 
since  cos  (0)  is  always  positive,  the  sign  value  of  a2  al¬ 
ways  reflects  the  sign  value  of  sin  (i|/): 

a2  =  cos  0  sin  i|i  (8.15) 

Therefore: 

$  =  acos  (^2e) '  (sisn  M)  (8-16) 

The  roll  ($)  angle  is  similarly  obtained: 

4>  =  acos  (^g)  ‘  (sign  [*s])  (8.17) 

9  Advantages  and  Disadvantages 

Quaternions  and  Euler  angles  have  their  own  ad¬ 
vantages  and  disadvantages.  Euler  angles  use  only  three 
components  instead  of  four  to  represent  orientation.  If 
one  were  to  send  quaternions  over  a  network  in  place  of 
Euler  angles,  as  has  been  proposed  for  SIMNET  (Burch- 
fiel,  1990),  network  traffic  would  increase.  However,  in 
cases  where  angular  rates  remain  constant  for  long  peri¬ 
ods  of  time,  by  extrapolating  orientation  updates  with 
quaternions,  it  would  be  necessary  to  send  an  update 
only  when  an  angle  rate  change  occurs. 


As  pointed  out  above,  quaternions  can  be  computed 
directly  from  the  dynamic  equations,  bypassing  the  com¬ 
putation  of  transcendental  functions  necessary  in  com¬ 
puting  Euler  angles.  If  each  transcendental  function 
costs  approximately  20  arithmetic  operations  then  the 
net  cost  for  deriving  a  rotation  update  using  the  Euler 
method  is  94  operations  (Burchfiel,  1990).  This  can  be 
compared  to  42  operations  using  the  quaternion 
method.  In  flight  simulation,  Euler  angles  are  sometimes 
necessary  for  use  in  other  simulator  functions  such  as 
cockpit  displays,  etc.  This  means  that  approximately  an¬ 
other  64  operations  are  necessary,  making  the  Euler 
method  slightly  more  computationally  efficient.  How¬ 
ever,  in  many  cases  it  is  not  necessary  to  compute  the 
angles  at  each  orientation  update.  Network  updates 
should  only  occur  when  a  rate  change  occurs.  Addition¬ 
ally,  radars,  gyros,  and  attitude  indicators  can  be  up¬ 
dated  at  slower  rates.  The  bottom  line  is  that  orienta¬ 
tions  can  be  achieved  very  efficiently  utilizing 
quaternions,  without  calculating  Euler  angles.  When 
Euler  angles  are  needed  for  other  aircraft  functions,  then 
quaternions  become  less  efficient,  depending  on  how 
often  they  need  to  be  computed  (Table  1). 

The  most  significant  advantage  of  quaternions  is  that 
no  singularity  exists  when  the  elevation  angle  (0)  passes 
through  ±tt/2.  In  the  Euler  method,  P  and  R  both  be¬ 
come  undefined  in  this  situation  due  to  division  by  zero. 
Techniques  exist,  however,  for  working  around  this  sin¬ 
gularity.  Truncating  values  as  ±tt/2  is  approached  will 
avoid  this  problem.  If  the  elevation  angle  is  truncated  at 
values  of  89.99  and  90.01  then  a  0.02  degree  rotation 
skip  results  (Goldiez  &  Lin,  1991).  Depending  on  the 
speed  of  the  program  and  the  rotation  rates  desired,  this 
may  not  be  noticeable.  However,  in  higher  fidelity  simu¬ 
lations,  where  a  slow  vertical  maneuver  is  executed,  it  is 
a  factor.  Regardless  of  the  significance  of  these  effects, 
however,  this  approach  has  the  disadvantage  of  intro¬ 
ducing  nonunique  values  for  Euler  angles. 

Numerous  aircraft  operate  autonomously  within  the 
prototype  simulation,  changing  very  little  in  angular  ve¬ 
locity.  When  this  simulator  is  eventually  networked  to 
other  workstations,  quaternions  will  provide  a  way  of 
forward  interpolating  rotations,  thereby  eliminating  the 
need  for  the  continued  transmission  of  update  packets.  If 


Table  I .  Efficiency  Comparison  of  Euler  and  Quaternion 
Methods 


Operation 

Euler 

Quaternion 

Derivation 

96 

42 

Creating  rotation  ma¬ 
trix 

20 

32 

Total 

116 

74 

Euler  angle  conver¬ 
sion 

0 

64 

Total  calculations 

116 

138 

updates  are  eventually  needed,  the  quaternion  rate  can 
quickly  and  easily  be  converted  into  Euler  angles. 

As  the  number  of  aircraft  increase  in  the  simulated 
world,  the  number  of  calculations  necessary  for  orienta¬ 
tion  updates  begins  to  multiply.  Since  only  the  currently 
piloted  aircraft  makes  use  of  the  Euler  angles  for  addi¬ 
tional  simulation  functions,  calculating  Euler  angles  is 
not  necessary  for  all  other  aircraft  in  the  simulation. 
Therefore,  it  becomes  more  efficient  to  utilize  quaterni¬ 
ons  for  defining  these  rotations,  saving  approximately  42 
arithmetic  operations  per  update  per  aircraft. 

1 0  Overall  System  Layout 

To  satisfy  the  prototype’s  basic  requirements,  the 
overall  structure  of  the  system  is  designed  as  shown  in 
Figure  13.  An  aircraft  data  file  makes  it  simple  to  create 
new  aircraft,  position  these  aircraft,  and  designate  their 
handling  characteristics.  It  is  also  used  to  initialize  the 
flight  parameters  of  the  autonomous  aircraft.  The  pro¬ 
gram  data  structure  contains  information  describing  the 
current  state  of  the  aircraft  and  its  design  specifications. 
The  aerodynamic  model  is  used  for  updating  the  piloted 
aircraft’s  body  rates.  The  nondynamic  model  also  out¬ 
puts  a  set  of  velocities,  but  unlike  the  dynamic  model,  it 
determines  these  velocities  in  accordance  with  a  prede¬ 
termined  script.  The  orientation  model  converts  body 
rates  to  position  and  orientation  in  world  space.  Euler 
angles  are  then  determined  for  the  piloted  aircraft  and 
used  to  update  cockpit  displays.  Autonomous  aircraft  do 
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Figure  1 3.  Prototype  flight  simulator  basic  structure. 


Flight  Record 
1 
1 

200.0 

-950.0 

300.0 

0.0 

0.0 

0.0 

3.0 

1 

Specification  Record 
4 
1 

27.5 
260.0 
10.8 
546.0 
8090.0 
25900.0 
29200.0 
1300.0 
0000.0 
8000.0 
0.03  0.3 

0.28  3.45  0.0  0.72  0.36 
0.0  -3.6  -0.38-1.1  -0.5 
-0.98  0.17 

-0.12-0.26  0.14  0.08  -0.105 
0.25  0.022  -0.35  0.06  0.032 
0.2618-0.5236  0.5236 


Ad  number/ 

/type  aircraft/ 

/airspeed/ 

/posx/ 

/altitude/ 

/posz/ 

/heading/ 

Anitial  angle  of  bank/ 

Anitial  gforce/ 

/status:0-ptioted,l-levelturn,2-g  controlled  turn/ 

/type  aircraft--  A4/ 

/jet  or  prop  jet:  1  prop:0/ 

AV 

/SI 

/c/ 

An/ 

Ax/ 

Ay/ 

Az/ 

Axz/ 

/max  thrust/ 

/mil  thrust  or  horsepower/ 

/CDo  /CDa 

/CLo  /CLa  /CLq  /CLda  /CLde 
/CMo  /CMq  /CMa  /CMda  /CMde 
/CYb/CYdr 

/CLb  /CLp  /CLr  /CLda  /CLdr 

/CNb  /CNp  /CNr  /CNda  /CNdr 

/deflection  limits  of  rud,  ail,  elevator  (radians)/ 


Figure  14.  Example  aircraft  data  file. 


not  require  the  calculation  of  Euler  angles  since  they 
need  not  display  instrument  readings  to  a  human  pilot 
and  therefore  bypass  this  function. 

1 1  Implementation  Details 

Data  records  within  the  aircraft  data  file  are  di¬ 
vided  into  two  categories,  flight  records  and  aircraft 
specification  records.  The  flight  record  contains  informa¬ 


tion  describing  the  position  of  an  aircraft  and  its  initial 
flight  parameters  (Fig.  14).  The  specification  record  con¬ 
tains  the  dimensional  characteristics  and  stability  coeffi¬ 
cients  describing  a  particular  aircraft.  By  manipulating 
the  data  in  the  specification  record,  one  can  change  an 
aircraft’s  basic  design  as  desired.  To  link  an  aircraft  to  a 
particular  set  of  specification  data,  all  that  is  necessary  is 
to  match  the  “type  aircraft”  information  within  the  two 
records.  This  system  allows  one  specification  record  to 
be  used  for  an  unlimited  number  of  flight  records. 
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int 

ID; 

—number  assignment 

int 

type; 

—aircraft  type 

int 

status; 

— piloted:0  or  autonomous  level  tum:l  climbing  tum:2 

typeptr 

Tptr; 

—pointer  to  aircraft  specification  data  structure 

float 

Forces  [3]; 

—forces  in  X,Y,Z  dir 

float 

Torques  [3]; 

-torques  around  X,YfZ  axis 

float 

linear_vel[3] ; 

-velocity  in  X,Y,Z  direction 

float 

angular_vel[3]; 

-angular  velocity  around  X,  Y  fZ  axes 

float 

linear_accel[3]; 

—linear  accelerations 

float 

angular_accel[3]; 

—angular  accelerations 

float 

sideslip; 

—sideslip  or  beta  angle 

float 

ang_atk; 

—angle  of  attack 

float 

d  ang_atk; 

—angle  of  attack  rate 

float 

lift; 

—total  lift 

float 

drag; 

—total  drag 

double 

Q[4]; 

—quaternion 

Matrix 

H; 

—direction  cosine  matrix 

float 

euler_angles[3]; 

— euler  angles  angles  in  radians  -  yaw,pitch,roll 

float 

pos[3]; 

—world  position  in  X,  Y,  Z 

float 

ref[3]; 

-look  direction 

float 

vel_world[3]; 

—velocities  in  world  position  -  X,Y,Z 

float 

gfor; 

-amount  of  g  force 

float 

rpm; 

—engine  rpm 

float 

elev; 

—flight  control  positions 

float 

eltrim; 

—elevator  trim 

float 

rud; 

-rudder  position 

float 

ail; 

-aileron  position 

float 

thro; 

—throttle  position 

int 

flaps; 

—flap  position 

int 

gear; 

—landing  gear  position 

Figure  1 5.  Aircraft  flight  data  structure. 

A  flight  data  structure  provides  a  global  source  of  in¬ 
formation  on  the  state  of  each  aircraft  in  the  simulation. 
The  information  contained  here  is  necessary  for  opera¬ 
tion  of  the  aerodynamic  and  orientation  models,  and  for 
updating  cockpit  displays  (Fig.  15).  Note  that  the  fourth 
item  in  this  structure  is  a  pointer  to  another  data  struc¬ 
ture  containing  aircraft  specification  data.  Not  shown, 
this  structure  contains  information  identical  to  that 
found  in  the  specification  record  of  the  data  input  file 
(Fig.  14).  Maintaining  flight  information  in  two  sepa¬ 
rate  files  saves  some  storage  space  by  allowing  more  than 
one  aircraft  to  use  a  single  set  of  specification  data. 

The  throttle  in  an  aircraft  is  the  pilot’s  primary  means 
to  control  the  engine.  As  such,  a  mapping  of  throttle 
position  to  engine  rpm  must  be  devised  that  incorpo¬ 
rates  delays  associated  with  engine  spool-up  characteris¬ 
tics.  In  large,  high-speed  simulators,  rpm  and  engine 
data  are  retrieved  from  engine-specific  lookup  tables. 
Because  tables  such  as  these  are  engine  specific,  the  fol¬ 
lowing  simpler,  generic  method  was  devised. 


Arpm  =  (rpmdesired  -  rpmcurrcnt)<p*  (11.1) 

where 

rpmdesired  =  throttle  position  x  throttle  gain 
dt  =  delta  time 

cp  =  engine  spool-up  gain  factor 
(inverse  of  time-constant) 

Since  applications  using  this  simulator  include  both 
jet  and  propeller  aircraft,  it  was  decided  that  allowances 
should  be  made  for  the  differing  characteristics  of  the 
two  types  of  engines.  Jet  engines  are  generally  rated  in 
terms  of  thrust  (lb),  while  propellers  are  rated  in  terms 
of  horsepower  (ft-lb/s)  (Anderson,  1989).  Since  the 
aerodynamic  model  uses  thrust  in  terms  of  lb,  it  can  use 
the  data  for  jets  direcdy.  However,  propeller  driven  air¬ 
craft  require  the  following  conversion: 
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aerodynamic  model 

read  aircraft  velocity 

time  step  for  aero  model  =  timestep_factor  *  velocity 
for  computed_time_step  loop 
do  aero  calculations 
update  aircraft  state  variables 

end  loop 

exit  aerodynamic  model 

Figure  1 6.  Algorithm  for  computing  time  step. 


T  = 


550r| HP  p 
Ft  p0 


(11.2) 


where 

t|  =  propeller  efficiency  (usually  around  .8) 

HP  =  engine  rated  horsepower 

p 

—  =  density  altitude  ratio  where  p0  is  the  density 

Po  .  , 

at  sea  level 

A  normal  aircraft  control  stick  exhibits  two  degrees  of 
freedom,  left-right  for  aileron  control  and  back-forward 
for  elevator  control.  Therefore,  control  inputs  from  the 
spaceball  were  limited  to  these  directions.  The  maximum 
deflection  of  the  control  stick  is  information  entered  via 
the  specification  records.  It  is  a  simple  procedure  to  read 
deflection  data  from  the  spaceball  and  linearly  map  it  to 
a  control  deflection  somewhere  between  ±max  obtain¬ 
able  deflection.  Rudder  deflection  was  not  simulated 
since  rudder  control  is  not  normally  used  in  jet  aircraft 
after  takeoff. 


1 2  Speed  of  Aerodynamic  Model 

The  aerodynamic  model  exhibits  a  tendency  to 
“blow  up”  if  the  time  step  between  aerodynamic  calcula¬ 
tions  becomes  too  great.  This  becomes  very  evident 
when  the  aircraft  speed  and  rotation  rates  increase.  The 
solution  to  this  problem  is  to  run  the  aerodynamic 
model  at  a  faster  rate  than  the  rest  of  the  system.  The 
trick  is  to  determine  how  fast.  “Smart”  integration 


schemes  can  be  found  in  the  literature  that  increase  or 
decrease  the  number  of  time  steps  according  to  the 
amount  of  numerical  divergence  present  in  the  inte¬ 
grated  values  (Press  et  al.,  1990). 

A  more  simple  method  is  used  in  this  program  to 
solve  for  this  problem.  Because  there  exists  a  direct  rela¬ 
tionship  between  an  aircraft’s  speed  and  its  tendency  to 
produce  a  numerical  instability,  the  time  step  was  ad¬ 
justed  in  direct  relationship  with  the  aircraft’s  airspeed. 
Prior  to  updating  body  rates  in  the  aerodynamic  model, 
aircraft  airspeed  is  measured  and  the  number  of  time 
steps  is  calculated  (Fig.  16). 

Increasing  the  time  step  does  not  decrease  perfor¬ 
mance  of  the  overall  system.  The  speed  of  the  aerody¬ 
namic  model  in  the  current  implementation  ranges  from 
17.6  to  18.8  msec.  Figures  17  and  18  illustrate  the  com¬ 
plexity  of  the  graphics  portion  of  the  overall  simulation 
model.  Running  at  approximately  120  msec,  the  graph¬ 
ics  pipeline  remains  the  limiting  factor  in  this  system. 


I  3  Conclusions  and  Future  Work 

The  techniques  presented  in  this  paper  have 
proven  to  be  an  effective  method  for  implementing  a 
graphical  dynamic  flight  simulation  on  a  matrix-based 
graphics  computer  in  real-time.  Like  most  research  and 
academic  projects,  this  aircraft  simulator  is  structured  to 
allow  for  the  addition  of  more  detailed  functionality. 
Current  work  includes  the  integration  of  a  weapons  de¬ 
livery  system  and  avionics  suite.  The  orientation  model 
functions  developed  during  the  course  of  this  paper  have 
become  part  of  the  C  program  library  at  the  Naval  Post¬ 
graduate  School,  thereby  providing  an  alternate  and 
more  flexible  tool  for  manipulating  solid  objects  in  a 
graphical  environment.  Integration  of  this  orientation 
model  into  other  dynamic  simulation  systems  is  also  un¬ 
der  investigation. 


Figure  1 7.  Looking  down  towards  runway. 


Figure  18.  Closing  in. 
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